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Abstract 



Motivated by experimental and theoretical work on autonomous oscillations in yeast, 
we analyze ordinary differential equations models of large populations of cells with 
cell-cycle dependent feedback. We assume a particular type of feedback that we call 
Responsive/Signaling (RS), but do not specify a functional form of the feedback. We 
study the dynamics and emergent behaviour of solutions, particularly temporal cluster- 
^ ing and stability of clustered solutions. We establish the existence of certain periodic 

I— —I clustered solutions as well as "uniform" solutions and add to the evidence that cell-cycle 

dependent feedback robustly leads to cell-cycle clustering. We highlight the fundamen- 
' tal differences in dynamics between systems with negative and positive feedback. For 

QQ positive feedback systems the most important mechanism seems to be the stability of 

^-H individual isolated clusters. On the other hand we find that in negative feedback sys- 

00 tems, clusters must interact with each other to reinforce coherence. We conclude from 

various details of the mathematical analysis that negative feedback is most consistent 
^-H with observations in yeast experiments. 

^ j Keywords: Inhomogcncous Feedback, Autonomous Oscillations in Yeast, Cell Cycle 

AMS Subject Classification: 37N25, 34C25, 34F05, 92D25 

1 Introduction 



1.1 Background 

In this paper we consider simple dynamical models of the cell division cycle. Specifically, 
consider a culture of n cells, in which the progression of the i-th cell is governed by the 
equation: 

—^ = \^a{xi,x). (1.1) 
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where Xi is the position of the ceh within the cycle and x denotes the state of all the cells 
in the culture. We will describe the dependence of a on Xi and x below. 

Our primary motivation for this model is recent theoretical and experimental work on 
Yeast Autonomous Oscillations (YAO) ([IZ], [SS], 0, [53], [M]), the periodic oscillations of 
physiologically relevant variables that have been reported and studied for over 40 years [9] 
[IIlllHllMlETllMlllQlllIlliailelllglEZI (and many others) . Different types of YAO 

have been called metabolic [55], glycolytic [2] or respiratory [23] oscillations. The control 
of oscillation and the regulation of yeast metabolism has been an important theme in the 
chemical engineering literature devoted to the efficient management of bioprocesses [TJ \2T\ 
[26l [53l [62] . These phenomena are of basic biological interest because they expose questions 
regarding the coordination of the cell cycle and metabolism, and interconnectedness of 
various cellular and genetic processes [H [31]. A correlation between YAO and the bud 
index was noted as early as [32[ I37| . However, it seems that the link between YAO and 
the Cell Division Cycle (CDC) was obscured by the fact that the periods of YAO always 
shorter than the CDC times (computed from dilution rate) and a relationship between YAO 
and the CDC seems to have been largely ignored. However in |31j . |55j and elsewhere, the 
correlation between YAO and CDC was again noted in genetic expression data. 



Figure 1: Phases of the yeast cell cycle. The Gl phase begins following cell division. The 
beginning of the DNA synthesis phase, S, coincides with budding. G2 is a second "gap" 
phase. The M phase is characterized by narrowing or "necking" between the parent and 
daughter cell; it ends in cell division. The hypothesized R region is the later portion of 
Gl and the signaling region S is in the S phase. That is, a large subpopulation of cells in 
the replicative S phase may promote or inhibit progression of cells approaching the Gl-S 
boundary. 

In |3] and [17] the authors proposed cell cycle clustering as a possible explanation of the 
interaction between YAO and the CDC. Figure [T] roughly illustrates the arrangement of 
phases of the cell cycle of yeast. We hypothesized that subtle feedback effects on CDC pro- 
gression could cause populations of cells to segregate into approximately CDC-synchronized 
subpopulations. Experimental bud index data reported in ^47J supported this picture. In 
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[3] the authors studied a few simple forms of (1.1) with the hypothesis that cehs in one part 
of the CDC may influence other cells in different parts of the CDC in different ways through 
various diffusible chemical products. We hypothesized that a large subpopulation of cells 
in the critical S-phase might effect metabolism production and the metabolites may in turn 
inhibit or promote cell growth in the later part of the Gl phase, thus setting up a feedback 
mechanism in which YAO and CDC clustering are inextricably intertwined. We showed 
analytically and numerically that differential CDC feedback such as this can robustly cause 
CDC clustering in the models. By clustering we do not mean spatial clustering (cultures 
that exhibit YAO occur in well-mixed bioreactors) , but groups of cells that are traversing 
the CDC in near synchrony. 
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Figure 2: Experimental time series from a continuous culture of budding yeast. Dissolved 
O2 percentage (green), bud index percentage (blue) and cell density (red) are plotted versus 
elapsed time. The average cell cycle period as calculated from dilution rate was about 400 
minutes. The plot shows clearly that the bud index (percentage of cells with buds from 
microscopy) and cell density (by flow cytometry) are both synchronized with the oscillation 
in the level of dissolved 02- 

Guided by these mathematical results, we verified the existence of clusters in two types of 
oscillating yeast using both bud index and cell density data. Some of the measurements 
from those experiments are shown in Figure [2j First we note that the cell cycle period, as 
calculated by the dilution rate, is approximately 400 minutes, and two O2 oscillations occur 
during this period, suggesting that there may be two clusters. Next, analyzing the figure, 
we see that approximately half of the cells are budding at times t = 50 and t = 250, while 
at t = 170 less than 10% of the cells are budded. Each budding event is accompanied by a 
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decrease in density (no cells are dividing) and followed by a sharp increase in cell density 
as these budded cells proceed through division. Note that at t = 170 since less than 10% 
of the cells are budded, most of the cells must be in the Gl phase of the cycle. When the 
bud index hits its next maximum at t = 250, approximately half of the cells must have 
budded. The other half of the cells must at that time still be in the Gl phase. As the 
cells that are budded then divide and the cell density increases, the cells that remained 
in Gl must still be in Gl since the bud index is again low. When the first group of cells 
has divided, the second group has been in Gl for at least 200 minutes. The next rise in 
bud index then must be due to these cells, since they have had time to mature, while the 
first group of recently divided cells clearly has not had time to reach budding again. Thus, 
these experiments show conclusively the existence of two clusters and that CDC clustering 
coexists with YAO. 



1.2 Modeling of the cell cycle and feedback 

In standard modeling the cell volume Vi{t) is a proxy for position in the cell cycle. This has 
justification for yeast in that milestones in the cell cycle, such as the onset of budding, are 
closely associated with volume milestones and thought to be causally related. Measurements 
show that the growth of a single cell is roughly exponential, so a first order approximation 
is that the volume of the i-th cell satisfies a linear differential equation 

- = cv.. (1.2) 

A frequent assumption on the growth rate c is that it does not depend on Vi i.e. it is 
independent of the cell's current state within the cycle and on other cells; it depends 
instead on the nutrients available and other environmental factors. Applying a logarithmic 
change of variables the growth law becomes dxi/dt = c, and by further normalizing both 
the coordinate x and time, the cell cycle can be represented by the unit interval [0, 1], and 
the equation of motion becomes dxi/dt = 1. In this simple model each cell reaches division 
(cytokinesis) at 1 when it returns (perhaps with its descendant cell) to and begins the 
cycle again. Note that a change to normalized coordinates does not depend essentially on 



the form of (1.2); if there is no interaction between cells and cells never stop growing, then 
one can change variables to the form dxi/dt = 1, with Xi(t) £ [0, 1]. 

A much more general model (again using normalized coordinates) is ( |1.1[ ). We proposed to 



consider forms of (1.1) where the cells in one portion of the cell cycle, S for signaling, may 
influence the growth rate of cells in a preceding portion that we term R for responsive. For 
example the R region may reside in the later portion of the Gl phase and the signaling 
region S may be the biological S phase (see Figure [T]). This is philosophically justified by 
the fact that the S phase is the most critical part of the CDC and the link between YAO 
and CDC may function to protect the integrity of transcription [9j. It is also known that 
growing yeast store carbohydrates, then liquidate them in the late Gl phase [T7]. The 
actual positions of the signaling and responsive regions within the biological cell cycle play 
no role in the mathematical analysis. In the rest of the paper we will use S to denote the 
signaling region and on the few occasions when we refer to the yeast's S phase we will do 
so explicitly. 
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Mathematically, the interval [0, 1], with the endpoints identified is a circle. On this circle, 
we can specify a positive direction as associated with the increasing direction on [0,1]. 
Distance between points x and y on the circle using these coordinates is given by the 
minimum of \x — y\ and 1 — \x — y\. 

Definition 1.1 Consider n cells whose coordinates are given by xi G [0,1], 1 identified 



with 0, and governed by an equation of the form (1.1). When a cell reaches 1 it continues 



at 0. We call such a system a RS-feedback system if: 

(HI) R is an interval that directly precedes another interval S, i.e. the last endpoint of R 
is the first endpoint of S, 

(H2) a{xi,x) vanishes except when Xi G R and there are some xj in S, 

(H3) < Umin < 1 + a{xi, x) < v-max for all Xi and x, 

(H4) feedback is monotone, thus adding a cell to S will increase the value \a{xi,x)\ for 
Xi € R, and, 

(H5) a{xi,x) is a smooth function for Xi in the interior or R and each xj in the interior 
of S, j ^ i and the one sided derivatives exist at the boundaries of R and S. 

By positive feedback we mean a is positive for Xi £ R if there are one or more xj in S. 
We define negative feedback in the analogous way. 

For the sake of definiteness we will specify throughout the rest of the paper: 

£• = [0,5) and R=[r,l), < s < r < 1. 

The final endpoint of i? is 1, which corresponds to 0, the initial endpoint of S. See Figure [3] 
below. 

Note that our restriction that R precedes S, while motivated by biological considerations, is 
not the only possibility. It is worth noting that in the reverse case, when R follows S, many 
of our results hold with the roles of positive and negative feedback reversed. In particular. 



this is true for Propositions 3.1, 3.2, 3.4, 4.1 4.2 and Corollary 4.3 



In [3j we considered two idealized forms of (1.1) with threshold-triggered feedback and 



found that with either positive or negative linear feedback, we robustly observe clustering. 
We also observed that the number of clusters formed is less dependent on the form of the 
feedback than on the sizes of the signaling and responsive regions, S and R. One aim of the 
present paper is to demonstrate that clustering typically occurs in a broad class of models 



(1.1), beyond the idealized situations investigated in [3j. The other aim is to study key 
differences between positive and negative feedback in the models in order to understand 
if one or the other is more favorable for clustering and thus more likely to be the process 
behind experimental observations. 

Note in Definition 11.11 that the number of cells in the culture is fixed at n. Thus we will 
consider neither death, harvesting or proliferation. In the oscillation experiments we are 
modeling there are in fact proliferation and harvesting (and an insignificant rate of death). 
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but they approximately balance when averaged over a cell cycle. Thus the expected number 
of cells descended from a single cell at any future cycle is approximately 1. Further, in the 
model we are considering, there is no distinction between the two cells resulting from a 
division and to keep track of both trajectories would be redundant. 



The differential equation in the general model ( 1.1 ) with RS feedback may have discontinu- 
ities. Thus we need to consider uniqueness and global existence of solutions. First consider 
that the equations may be discontinuous only when a variable is at the boundaries of R and 
S i.e. at the hyper-surfaces given by xj = 0, s, r for some j. We obtain uniqueness since by 
(H3) any solution crosses a surface of discontinuity non-tangentially with non-zero speed. 
Global existence follows from boundedness of the vector field (H3). 

A version of the model that we will study in sections 6 and 7 is: 

dxj,_\l, Xi^R 
dt ~ \l + /(/), if Xi£R 

where 

^ #{J ■■ X, £ S} ^ #{j:x,e[0,s)} ^ 
n n ' 



i.e. I is the fraction of cells in the signaling region. The "response function" /(/) in (1.3) 
must satisfy /(O) =0 and be monotone, but can be non- linear, for instance sigmoidal (S- 
shaped). We will see that models of this form, while fairly general, can be studied in some 
detail. In [15] we studied some cases where / is a linear function. 

Understanding the CDC at the genetic and biochemical level is a topic of intense interest and 
progress has been made in identifying the agents and the nature of relationships between 
them Uni [TTl [Sni ED ES] . Our approach uses a "caricature" of the cell cycle, rather 
than detailed modeling, and this simplification demands justification. First, we wish to 
deal with individual cells in a population-wide phenomena. If details within each cell are 
considered, then the dimensions of the phase space would be extremely large and results 
would be difficult to obtain. Second, our understanding of the details of the cell cycle 
and its relations with other processes is not complete and even if the general nature of 
relationships were well-understood, the resulting set of differential equations would contain 
many parameters, e.g. rate constants, that could only be estimated. With our simplified 
model which is based on biological insight, we hope to obtain general principles that will 
inform further detailed investigations. 

The approach in part of this work is basically that of "phase oscillator" models, e.g. Ku- 
ramoto equations, in which details of each individual actor are projected onto a simple 
phase space and emergent population behaviours are studied. In fact if /(/) is linear then 
our model can be put into the Kuramoto form by integrating over the cell cycle for each 
pair of cells and adding the effects (see |33] p. 65-67). This derivation fails in general for 



RS models (1.1) or (1.3) since the effects of cells in general are not additive. 

We note here another modeling simplification that we have implicitly made; namely one 
might more accurately model the feedback term as a[xi,z) where z \s a, vector variable 
representing all substrate factors that contribute to growth rate and z itself is coupled with 
x [25]. Dropping the z variable can be justified if the time-scale of the dynamics of this 
variable is significantly shorter than the time-scale of the CDC. 
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2 Clusters, Gaps and Isolation 



In this section we begin to study the existence and stabihty of periodic "clustered" solutions 



for both positive and negative feedback systems (1.1) and in later sections we point out 
crucial differences between these two types of feedback. Reducing to the study of clustered 
solutions is of practical interest since clusters appear in experiments with YAO. It also 
limits the dimensions of the problem to a manageable size. This strategy has proven 
indispensable in many fields; for instance in fluid dynamics, insight is obtained by studying 
finite dimensional vortex equations rather than the full Navier-Stokes partial differential 
equations 



Definition 2.1 By a cluster, we will formally mean a group of cells that are completely 
synchronized in the CDC. 

Note that RS-feedback systems as defined have the symmetry of globally coupled networks 
with identical nodes; namely, the vector field is equivariant with respect to the group of 
permutations of coordinates. This symmetry implies that any cells that initially share the 
same phase keep the same phase as time evolves. The simplest trajectory consists in taking 
all cells initially in the same phase. We have a single cluster C (synchrony) that generates 
a periodic solution that runs at velocity 1 around the circle {i.e. Xi{t) = Xj(0) + t mod 1 
for every i and all t > 0). 

Definition 2.2 By the gap between two clusters or cells at Xj_i and xi we mean the open 
interval on the circle from Xi-i to Xi, in the direction of the flow that contains no cells and 
has width Wi = Xi — (This can be made precise using a lift to the real line.) 

Note that if there are only two clusters xq and xi in the system, then there are two gaps 
between them and each gap has a width. 

It follows from our assumed coordinates that if two clusters are in S and R then the distance 
between them on the circle is less than + \S\ = 1 — r + s. Here \R\ denotes the length 
of the interval R = [r, 1), which is 1 — r and IS*! = s denotes the length of the interval 
S = [0, s). See Figure [sj We say that a cluster of cells is isolated if there are gaps between 
the cluster and any other cells on either side of length at least + IS"! and strictly isolated 
if the widths of gaps are more than + \S\. This terminology is motivated by the fact 
that strictly isolated clusters cannot exert feedback on cells outside the cluster, or have 
feedback exerted upon them from outside. While we consider only clustered solutions in 
the strictest sense, in real cultures individual cell differences will lead to a weaker form of 
clustering. For clarity we will refer to such a weakly clustered subset of cells as a group of 
cells. In Figure [3] the solution in the simulation has formed groups. 

The following definition will play a large role in the analysis of the model. 
Definition 2.3 Define 

M^Ld^i + i^ir'j- (2.1) 

M is the maximum number of isolated clusters that can simultaneously exist, given the sizes 
\R\ and \ S\ of R and S. 
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Figure 3: Our coordinate representation of the cell cycle and (weakly clustered) groups of 
cells from a negative feedback simulation with n = 200 and parameter values s = .25 and 
r = .75. Positions of individual cells are denoted by red asterisks. In this coordinate system 
the S region is the interval [0, .25) and the R region is [.75, 1), where 1 is identified with 0. 

Here [xj denotes the floor function, that is, the greatest integer less than or equal to x (e.g. 
[2.1J = 2). In Figure g when {R\ + |5|)-^ is in the range [3,4), M is 3. 

Proposition 2.4 For any RS feedback model and any positive integer k < M there exist 
periodic solutions consisting of k isolated clusters that do not interact. 

Proof: For k < M consider the solution with initial conditions: xq = 0, xi = ^, 
rcfc = We claim that this is such a periodic solution. Notice from the definition of M 



Since the distance between any two consecutive clusters x and y is initially d{x,y) = 1/k, 
no two clusters can be in i? U 5 simultaneously. Thus no feedback will occur and thus the 
distance between clusters will not change. 

By the same reasoning any initial condition of k < AI clusters, where all pairs of clusters 
satisfy d{x,y) > \R\ + will also lead to a periodic solution where all clusters move 
indefinitely with speed 1. □ 

Conversely, if more than M clusters exist, then at least two of them are within a distance 
+ [S"! of each other, and while the first of these clusters lies in the signaling region, it 
will exert feedback on the second cluster for a non-empty interval of time. 



that 



M<{\R\ + \S\) 



-1 



and so 
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Note that a solution consisting of strictly isolated clusters can be at most neutrally stable 
(not asymptotically stable) since moving a cluster to the left or right still produces an 
isolated cluster. 

In Figure |4] we plot the results of numerical simulations which compare the number of 
clusters that formed with the maximum number of possible isolated clusters. The number 



of cells n was 5000 and the model was taken to be that in equation (1.3) with / linear. 
Specifically, for the positive feedback simulations /(/) = .61, and /(/) = —.6/ in the 
negative feedback simulations. There was also stochastic noise added to the equations 
with level a = 10~^ in order to demonstrate robustness. In the plot the x-axis represents 
{\R\ + 15*1)"^ which was varied from 1 to 6 in one hundred increments. In the simulations 
\R\ and \S\ were taken to be equal. The system was integrated for up to 100 cell cycles 
(with 50 steps per cycle) to check for the formation of clusters (but clusters usually formed 
long before the 100th cycle). To test for clustering, we produced histograms with 120 bins 
and visually inspected these for clustering. If clusters did not clearly form for a given value 
of (|i?| + IS"!)""*^, then that data point is plotted on the graph as a zero value for A^. 
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Figure 4: The number of clusters that form in simulations compared with M 
\S\)^^\ , the maximum number of isolated clusters given R and S. 



\R\ + 



The most striking feature of the plots in Figure |4] is that for positive feedback the number 
of clusters formed is always less than or equal to M, but for negative feedback the number 
of clusters formed is always greater than M. Notice also that positive feedback always 
produced clustering, but this was not the case for negative feedback. For negative feedback 
there are parameter values where no clusters form. Finally, it is worthy of note that for 
negative feedback, there are no occurrences of one-cluster solutions. Analysis in the next 
sections will shed light on these observations. 



9 



Besides solutions consisting of k < M isolated clusters, RS-feedback systems have other 
periodic solutions. One of these consists of all n cells spread along the cycle as uniformly as 
possible. We will define a uniform solution to be a trajectory for which the coordinates 
satisfy the following relation for some time d > 0: 

= Xj-i-i(O) for all z = 0, ...,n — 2, and = xo(0) modi. (2-2) 

Since the velocity of cells in the complement of R is precisely 1, it follows that for such 
a solution the cells in R'^ will be uniformly distributed with inter-cell distance d. In the 
process of the construction of uniform solutions we also produce many other periodic k 
cluster solutions. 

Suppose that k divides n and and k sets of n/k cells are initially synchronized. Then we 
may greatly reduce the dimensions of the differential equations by considering only the 
positions of the k clusters which we may denote by {xQ^t), xi{t), . . . , In particular 

we will prove that there always exists a solution of k clusters that satisfy: 

= Xj+i(0) for all i = 0,...,k — 2, and Xk-i{d) = xo{0) modi. (2-3) 

We will refer to such solutions as cyclic k cluster solutions. 

Proposition 2.5 There exists a uniform solution of any RS-feedback system. If k is a 
divisor of n, then a cyclic k cluster solution exists consisting of n/k cells in each cluster. 

We defer the proof of this proposition until Section 5.2. We observe that we do not have 
uniqueness of the uniform and cyclic clustered solutions. 

For n large, as in the application in mind, then we expect k cluster solutions for k « n 
to exist even if k does not exactly divide n. For example hyperbolicity (linear stability or 
instability) of the k cyclic solutions when k divides n would imply that such solutions exist 
for all n' n. 

The uniform solution has an analogue in PDE models of the cell cycle which we note in 
Appendix A. 



3 Positive Feedback Systems 

In a RS model with positive feedback, first notice that a group of cells that is isolated will 
remain isolated. Further, positive feedback has a focusing effect on an isolated group. 



Proposition 3.1 In a general RS model (1.1) with positive feedback, suppose that a so- 



lution has a gap between two adjacent cells Xi-i and Xi of width greater than or equal to 
\R\ + \S\. Then the width of this gap will never decrease. In particular, an isolated group 
or isolated cluster will remain isolated indefinitely. 

If there are only two clusters in the system, then this proposition applies to either of the 
two gaps that has width at least \R\ + l-SI. 
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Proof. Suppose that two consecutive cells and Xi are separated by a gap of width 
it'i > + l-S*]. Since feedback is positive, the cell Xi always moves at a speed of at least 



1. The cell speed of being governed by (1.1), will be exactly 1 whenever Wi is greater 



than or equal to + IS*!, since there will be no cells in S when Xi-i is in R. We then find 
that the time derivative of Wi is non-negative when Wi > + 15*1. Therefore, the gap can 
never decrease. It follows immediately that a group or cluster of cells that is isolated will 
remain isolated. □ 



Proposition 3.2 In a general RS model (1.1) with positive feedback, suppose that a group 



of cells with width w less than \R\ + \S\ is isolated. Then the width of the group will converge 
to zero as t ^ +oo. 



Proof. By the previous proposition the group will remain isolated for all future time and 
we may consider this group of cells as a decoupled sub-system. Without loss of generality we 
may renumber the cells in the group so that they have coordinates: xi, X2, • • • ,X£, (ordered 
in the direction of the flow). For i = 1,. . . ,i — 1, denote by Wi = Xi — Xi the width of 
the interval from Xj to X£. By assumption w^-i < Wi < wi < \R\ + \S\. Observe that this 
condition ensures that each Xj, i = 1, . . . — l will experience some acceleration every time 
it passes through R, unless Xi = X£. This implies that if uJi{t) is non-zero, then Wi{t) will 
decrease each time the group passes through R and 5. 

Since the group will remain isolated, the cell Xi will always move with speed 1 and Wi{t) will 
never increase. Since each Wi{t) is non-increasing and bounded below by 0, it must have a 
limit . Now consider a solution with an initial condition such that Xj(0) = x^(0) — w°° 
for 1 < i < i. By a standard argument using continuous dependence of solutions on initial 
conditions, for this solution each Wi{t) will be identically w°°. This implies that wf° must 
be zero, since, from above, a non-zero Wi{t) < -|- IS*! must decrease during each cycle. □ 

In the next result we discuss stability, for which we need the concept of distance and 
neighborhoods in phase space, which for the models we are considering is the n-torus, T", 
where n is the number of cells. Note that on T" there is a natural metric (distance) defined 
by the maximum of the (mod 1) coordinate differences. 



Proposition 3.3 In a RS model 



1.1) with positive feedback, the set of strictly isolated 
A solution consisting of k > 2 strictly 



cluster solutions is locally asymptotically stable, 
isolated clusters is neutrally stable (stable, but not asymptotically stable) within the set of 
solutions with k clusters. 



Proof. First observe that an e-neighborhood of a configuration consisting of isolated 
clusters consists of groups of cells within e of the original clusters. If the original clusters 
are strictly isolated, then we may make e small enough that the groups are also strictly 
isolated. By Propositions 3.1 and 3.2 each of these groups will remain isolated and converge 
to a cluster. Thus a solution starting at any initial condition within a neighborhood of the 
set of strictly isolated clusters will asymptotically approach the set. 

The second part of the claim follows since, if a strictly isolated cluster is moved a small 
distance, then it is still strictly isolated. Thus a small perturbation of a solution consisting 
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of k strictly isolated clusters also will be a solution consisting of k strictly isolated clusters. 
The distance between the two solutions will remain constant for all future time and thus 
they are stable, but not asymptotically stable. □ 

Points near the set of isolated cluster solutions will converge to the set, but individual 
solutions are only neutrally stable with respect to perturbations inside the set. 

Since the set of clustered solutions is locally stable, it must have a basin of attraction and 
it is natural to ask how big the basin is. In simulations for positive feedback systems, the 
basin seems to include almost all initial conditions. In the next proposition we see that the 
basin of attraction extends far beyond a small neighborhood of the set. 



Proposition 3.4 Suppose that a solution x{t) in a RS model (1.1) with positive feedback 
has at least one gap of width greater than or equal to \R\ + \S\. Then the solution will 
converge to a periodic solution consisting entirely of isolated clusters. 



Proof. We will call a gap large if its width is greater than or equal to \R\ + At time 
to, the cells may be grouped into a minimum number of groups in which there are no large 
internal gaps. The number of such groups is the same as the number of large gaps. Note 
that each such group is isolated, and the number of such groups cannot be larger than 
M = [{\R\ + |5|)~^J . Consider one such group. Since it contains no large internal gaps and 
it is isolated, during passage through R the last cell must be accelerated by the presence 
of at least one cell in S and so its speed is sometimes greater than 1. On the other hand, 
since the group is isolated it will remain isolated by Proposition 3.1, and the lead cell will 
travel indefinitely with speed 1 by the same argument that appears in the proof of that 
proposition. 

Thus during one passage through R the distance between the lead cell and the final cell 
in the group must decrease. If this group continues to have no large internal gaps, then it 
follows that the width of the group will continue to decrease. By an elementary argument 
(as for Proposition 3.2), the width will converge to zero; in other words the group will 
converge to a cluster. Otherwise, if a large internal gap develops then the cells that are 
separated by the gap will be isolated from each other and thus form two isolated groups. 
When this occurs the number of isolated groups will increase. Since the number of isolated 
groups is bounded above, large internal gaps may form only a finite number of times and 
thus eventually we have a fixed number of groups that never develop large internal gaps 
and each of these converges to an isolated cluster. □ 



4 Negative Feedback Systems 

The key observation is that for negative feedback, isolated clusters are not stable. This 
is because as a group of cells crosses the R-S boundary all cells of the group are delayed 
except the lead cell, which moves with unit velocity, causing the group to spread. 



Proposition 4.1 In a RS model (1.1) with negative feedback, a solution consisting of 



strictly isolated clusters is locally unstable. 
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Proof. Denote by x*{t) a solution consisting of strictly isolated clusters. First observe 
that under the condition of strict isolation the gaps between clusters are all larger than 
|-R| + |5|, and so any sufficiently small perturbation of the clusters consists of groups that 
are still isolated. In any neighborhood of any configuration with an isolated cluster, there 
is a configuration where Xj 7^ x*{0) for any i. Note in fact that a local coordinate system in 
a neighborhood of x*(0) is {xi — x*}, < i < n — 1. Now let x{t) be a solution with initial 
condition x{0) that differs from x*{0) in only the i-th coordinate. When this cluster passes 
through the boundary from R to S the separation between Xi and the rest of the cluster to 
which it belongs will increase. If we further let the perturbation be sufficiently small, then 
the cluster will remain isolated. 

Now recall the definition of stability: given any e > 0, there exists 6 > such that any 
solution x{t) starting within a S neighborhood of x*{t) will remain indefinitely within an e 
neighborhood of x*{t). Let eo be the largest e so that any x(0) within an e neighborhood of 
x*{0) will consist of isolated groups and let e = eo/2. If x(0) is as in the previous paragraph 
and is arbitrarily close to x*(0) then the distance between x{t) and x*{t) will continue to 
increase on each unit time interval as long as x{t) continues to consist of isolated groups. 
Therefore it follows that x{t) will eventually be outside of an eo/2 neighborhood of x*(t). 
Thus x*{t) is not stable. □ 

Note that we have not proved linear instability (when the derivative of the mapping at the 
periodic point is greater than one in absolute value which implies that nearby orbit are 



repelled exponentially). Linear instability can be shown for systems of the form (1.3). 

It follows that in order for clusters to remain coherent under small perturbations in negative 
feedback, they must not be isolated, i.e. the gaps between them must be less than + 
and so the number of stable clusters must be at least M + 1. This is clearly confirmed in 
simulations. Further it seems that sometimes M + 2 is a stable number of clusters. We will 
show in §7 that a A; = M + 1 cluster cyclic solution is stable for negative feedback of the 



form (1.3) and for some open sets of parameter values. 

In the following proposition we see that interacting clusters tend to spread out from each 
other as far as possible. 



Proposition 4.2 In a RS model (1.1) with negative feedback suppose that two clusters are 
within \R\ + \S\ of each other, hut are isolated from other cells (non-empty). If they remain 
isolated from other cells, then the gap between the two clusters will increase and converge 
to \R\ + IS"!. In the case that the two clusters contain all the cells in the system, if one 
gap has width less than \R\ + \S\ and the second gap has width > \R\ + \S\, then as long 
as the width of the second gap remains greater than + |5| the first gap will increase and 
approach + |5| . 



Proof. If the gap width (or smaller gap width in the case of only two clusters) is less than 
\R\ + IS"! then each time the second cluster passes through R, the first cluster will be in S 
for a non-empty interval of time. During this time interval, the cluster in R will experience 
deceleration, and the width of the gap will increase during the passage through R. If the 
cluster pair remains isolated, then the distance will be preserved through the rest of the 
cycle. Thus the distance between clusters will increase during each cycle. The distance is 
bounded above by + IS*! and so, by a standard argument, the sequence of distances thus 
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generated will converge to + IS"]. □ 



Corollary 4.3 Suppose that there are k < M clusters in a RS system (1.1) with negative 
feedback. Then the solution will converge to a periodic clustered solution with isolated 
clusters. 



5 Dynamics of clusters via return maps 
5.1 The return map for a clustered system 



Let us continue to consider the most general model ( 1.1 ) with RS feedback. Let a population 



of n cells be organized into k equal clusters and let the clusters be labeled by a discrete 
index i G {0, . . . , A; — 1} so that x = {xi\^~Q represent clusters of n/k cells each. One can 
assume that all coordinates Xi{{)) of the k clusters are initially well-ordered as 

= xo(0) < . . . < Xi(0) < . . . < 1, i = 1, . . . , A: - 1. 

This ordering is preserved under the dynamics (this can be well-defined using the orientation 
of the circle and xo(t) as a moving reference point). Moreover, the first coordinate xq must 
eventually reach 1, i.e. there exists such that xo{t-^{) = 1. Thus the set xq = defines a 
Poincare section for the dynamics and the mapping 

(xi(0), X2(0), . . . , Xfc-i(O)) ^ (xi(tR), X2(tR), . . . , a;fe_i(tR)) 

defines the corresponding return map. 

Starting from t = 0, compute the time ti that Xfc_i needs to reach 1 and compute the 
location of the remaining clusters at this time. Define a map F by 

F : (xi(0), X2(0), . . . , Xfc_i(0)) ^ (xo(ti), xi(ti), . . . , Xfc_2(ti)). 

Notice that xo(ti) = ti by assumption on xo(0). An illustration of F in the case A; = 3 is 
given in Figure [5} 

1^0 1 
(j) o — o 1 I — y I o o — cj) 

Xo X\ X2 X'l x'2 x'q 

Figure 5: An illustration of the map F with fc = 3; F(xi,X2) = (x'^jX^. 

Now the time ti +t2 that Xk-2 needs to reach 1, together with the population configuration 
at t = ti + t2, follow by applying F to the configuration (xo(ti), xi(ti), . . . ,Xfc_2(ti)). By 
repeating the argument, the desired return time tR is given hy t^ = ti + t2 + • • • + tk and 
the desired return map is F^'. Therefore, to study the dynamics one only has to understand 
the first map F. 

Note the Kuramoto derivation (in which one averages over a cell cycle for each pair of cells 
then adds the effects ([33| p. 65-67)) is a specific way to calculate F'' that is applicable 
when the feedback is additive. 
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We will first consider general properties of the map F and then we will compute and analyze 
it for feedback of the form (1.3) in the simplest case, k = 2. We emphasize that the case 
k = 2 is also perhaps most important since it is the first to have been directly observed. 



5.2 General properties for arbitrary k 

We may regard F as a continuous piecewise affine map of the {k — l)-dimensional simplex 

< 2;i < X2 < • • • < x/c-i < 1 

into itself. (Although the boundaries and 1 are identified in the original flow, in the 
analysis here, we consider them as being distinct points for convenience.) 

On the edges of the simplex, F has relatively simple dynamics. Indeed, if initially all 
coordinates are equal, then they must all reach the boundary 1 simultaneously. In other 
words, on the diagonal (xj = x for all i), we have F{x, . . . ,x) = (ti, 1, . . . , 1) where ti 
depends on r,s and x (for x = 0, we have ti = 1 independently of r and s). Moreover, 
starting with x^-i = 1 implies ti = which yields 

F{xi, Xk-2, 1) = (0, XI,..., Xk~2) 

whatever the remaining coordinates xi, . . . ,Xk-2 are. As a consequence, the edge 

{(x,l,...,l) : XG [0,1]} 

is mapped onto 

{(0,x,l,...,l) : xe [0,1]} 

which is mapped onto {(0,0, x, 1, . . . , 1) : x G [0, 1]} and so on, until it reaches the edge 
(0, . . . ,0,x), which is mapped back onto the diagonal (after k iterations). 

A particular orbit on the edges is the /c-periodic orbit passing the vertices, and which 
corresponds to the single cluster of velocity 1 in the original flow, namely 

(0, . . . , 0) ^ (1, . . . , 1) ^ (0, 1, . . . , 1) ^ (0, 0, 1, . . . , 1) ^ . . . ^ (0, . . . , 0, 1) ^ (0, . . . , 0). 

Geometrically, the corners of the simplex are cyclically permuted by the map F. It follows 
that the k — 2 dimensional simplexes (faces) that make up the boundary of the k — 1 simplex 
are also cyclically permuted by F. This implies that F cannot have a fixed point on the 



boundary. Propositions 3.2 and 4.1 tell us this orbit (which represents a single isolated 
cluster) must be asymptotically stable for positive feedback and unstable for negative. 
Since there are two points of a periodic orbit at each boundary of every edge (which are 
themselves globally /c-periodic 1-dimensional sets), and since these points are either both 
stable or both unstable, there must be at least one other A;-periodic orbit on the edges with 
coordinates between and 1. Whether this orbit is unique might depend on parameters. 

Proof of Proposition 2.1. Since the simplex is a convex and compact invariant set under F, 
the Brouwer fixed point theorem implies the existence of a fixed point. Since the boundary 
cannot contain any fixed point, the fixed point is in the interior. Note that the fc-cluster 
cyclic solutions are fixed points of F and vice versa. The uniform solution follows by taking 
k = n. □ 
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6 Dynamics for k = 2 



In this section we study the dynamics of 2 cluster systems for the model ( 1.3 ). Studying the 



behavior in the cases of a small number of clusters is not just a matter of convenience, but 
is important from the perspective of applications since presumably only a small number 
of clusters can form (for fixed n, more clusters implies that each cluster contains fewer 
cells and thus can exert less influence) and be observable (smaller clusters would produce 
smaller oscillation in metabolites and other chemical agents). In the experiment reported 
in Figure [2] there are 2 clusters. 

6.1 The map F 



Consider (1.3) where / is a monotone function. In the case k = 2, since only one cluster 



can exert feedback on the other, RS-feedback (1.3) simplifies to: 



dxi 
~dt 



l + fil) if Xi G [r,l) and Xj G [0,s),j ^ i 
1 otherwise 



Let a = f{-^) for notational simplicity. 

For k = 2, F is defined on the interval [0, 1] and is determined by xo(ti) = ti where ti is 
the time at which xi{t) reaches 1. When regarded as a function of xi only, its explicit form 
depends on the parameters r and s. There are two cases depending on the relative sizes 
of the signaling and responsive regions, specifically on the size of (1 + a)s with respect to 
1 — r. We have put the details of the computation and analysis of F{xi) in Appendix B. 

In the case where r + (1 + a)s < 1, we obtain that F is a continuous decreasing map: 



F{xi) 



In the case where r 



1 



1 — Xl 
(1 + a)xi + a(r 



1 — Xl 



- as 

Xl) 



^(1 

{1 + a)s > 1, the map F is as follows 



if < Xl < r — s 

if r — s < Xl < r 

if r < Xl < 1 — (1 + as 

if 1 - (1 + a)s < Xl < 1. 



Fixi] 



1 — Xl 

1 — (1 + Q;)xi + a{r — s) 
r-xi + -^il-r) 

aTT(l-^i) 



< Xl < r 
r 



if 
if 

a+l ^ a+l' 

if r < Xl < 1. 



5 < Xl < ^ + 

s < Xl < r 



Calculating the full return map -F^ is prohibitively complicated, but in Appendix B we use 
these two possible forms of F to analyze the dynamics. In the next section we summarize 
the results. 



6.2 Analysis of the dynamics 



In Appendix B we find only four distinct types of dynamics; two for positive feedback and 
two for negative. 
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• Positive feedback: 

1. There is a unique unstable fixed point for F^. 

2. There is an interval of fixed points for F'^. 

All other orbits are asymptotic to the boundary (merger of the two clusters). 

• Negative feedback: 

1. There is a unique stable fixed point for F'^. 

2. There is an interval of fixed points for F^. 

All other orbits, except the boundary points, are asymptotic to the stable fixed point 
or the interval of fixed points. 

These possibilities, for some specific parameter values, are illustrated in Figure |6j 




Figure 6: Plots of the Poincare map F"^ in the case k = 2 for various parameter values 
{a = /(^)). Clockwise from top left: interval of fixed points under negative feedback, 
unstable fixed point under positive feedback, interval of fixed points under positive feedback, 
stable fixed point under negative feedback. 
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An interval of fixed points, we observe in the Appendix B, can occur not only because the 
two clusters may be isolated from each other, but also in certain other situations. Namely 
if either: 

• xq is in S for the entire time that xi is in R, or, 

• xi is in R for the entire time that xq is in S, 

then the unique fixed point of F is neutral and contained in an interval of neutral period 2 
points (fixed points of -F^). 

In [15] we present similar computations for a subset of parameter values with three clusters 
and positive linear feedback. There we compute F where F^ is the full return map. The 
results there are similar to those reported here; for all the cases examined the three cluster 
cyclic solution is either unstable or in the interior of a set of neutral periodic solutions 
(period 3 points of F) . No other periodic orbits were detected and all other initial conditions 
tend to two cluster or one cluster periodic solutions (on the boundary of the domain of the 
map F. 

A summary of our studies of A: = 2 and k = 3 cluster systems is as follows. If the system 
has positive feedback, then many initial conditions lead to a single cluster, but if the initial 
condition begins with 2 or 3 clusters, or close to such, then these clusters might persist 
depending on the parameters and initial conditions {e.g. if the clusters are isolated). If the 
system has negative feedback, then there may be solutions with 2 or 3 clusters (depending 
on the parameter values) that are stable within the set of clustered solutions. One cluster is 
never stable under negative feedback. Biologically, synchrony is likely to appear in systems 
with positive feedback and clustering in systems that have negative feedback. 



7 Cyclic M + 1 Cluster Solutions 



Again consider the model (1.3) of RS feedback. Recall that M = [{\R\ + IS"]) is the 



maximum number of clusters that can exist without mutual interactions. In this section we 
consider the cyclic solutions consisting of A; = M + 1 clusters, with coordinates xq, ...,Xk-i, 



in the dynamics corresponding to (1.3) 



We again denote the signaling region S by [0, s) and the responsive region i? by [r, 1). 



Proposition 7.1 Consider RS feedback of the form (1.3). For any < s < r < 1, there 
is a cyclic solution consisting of k = M + 1 equal clusters of the form xq = 0, xi = 
d, Xk-i = {k — l)d, for some d > 0. Denote f3 = /(^). // 



1 A +/3r 
^^k V 1 + /3 



and r > 



k-1 
k 



(1 - s/3) 



(7.1) 



then the fixed point is unstable for positive (3 and stable for negative (3. Otherwise, the 
solution is neutrally stable. Both stability results are within the set of k-cluster solutions. 



(Note that for k = 2 the first (resp. second) inequality in (7.1) corresponds to the re^ 



quirement for the second line segment in Figure [8p (resp. [8p) to meet the diagonal, as 
k = M + 1 = 2 forces r - s < |) 
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s 



Figure 7: Regions of parameter space for the k cyclic solutions with k = M + 1. Each 
diagonal band contains the parameters for a specific k and is partitioned into three cases. 
In case I the /c-cyclic solution is unstable for positive feedback and stable for negative. In 
cases II and III the A:-cyclic solution is neutral and is contained in a set of neutral period 
k solutions. In this plot /3 = /(^) is taken to be 1/k. 



Note that if the solution is neutrally stable, then the fixed point of F must be contained in 
the interior of a set of period k points that are also neutrally stable, since the return map 
is piece-wise affine. Note that this is consistent with the results for two and three cluster 
systems. 



The proof of Proposition 7.1 is contained in Appendix C. 



8 Discussion 



For both positive and negative feedback, the maximum number of non-interacting clusters, 
M, (see Definition 2.3), plays a large role in the dynamics of the system. 

Our main result for positive feedback is that the strictly isolated clustered solutions are 
neutrally stable (within the set of k clustered solutions, k < M), while the set of strictly 
isolated clustered solutions is locally asymptotically stable (in the full space). The basin 
of attraction of this set extends beyond a small neighborhood of the set. For (1.3) with 
positive feedback and some regions of parameters we have proved that the periodic solution 
of A; = M + 1 equal clusters is completely unstable. In simulations all initial conditions 
seem to lead to clustered solutions with k < M. Given that these are only neutrally stable 
within the set of clustered solutions, solutions may quickly approach clustered solutions, 
but small perturbations of the model may eventually cause the clusters to merge into one 



19 



cluster (full synchronization). 

For negative feedback we observed that isolated clusters are unstable. This implies that 
clusters need interaction with other clusters to remain coherent. When clusters form, the 



number of clusters k is at least M + 1. For the model (1.3) and under some restrictions on 
the parameters, we have proved that the k-cyclic periodic solution of A; = M + 1 clusters is 
stable. In numerical experiments, we observed that clusters sometimes do not form under 
negative feedback. This perhaps indicates that complicated bifurcations occur in between 
regions of clustering. 

From simulations in the realistic case n » M the uniform solution seems to be unstable for 
both positive and negative feedback. Analytic treatment of stability of uniform solutions 
is a challenge for future work. 

Perhaps the key observation from the mathematical models is that positive feedback typ- 
ically leads to synchrony, while negative feedback systems tend to clustering. In practical 
terms, if a positive feedback mechanism similar to ours were in effect then one would expect 
to encounter synchronization, at least some of the time. In yeast experiments a single cluster 
(fully synchronized behavior) seems to be impossible to sustain in the laboratory [6l [5^ [59] 
and in YAO experiments there is clearly no CDC synchrony. This creates a strong suspi- 
cion that the causal mechanism underlying clustering and YAO involves negative feedback 
rather than positive. 

Acetaldehyde [2J, ethanol, oxygen and other membrane-permeating metabolites [25j have 
been conjectured as possible mediators of the YAO. It is natural that these same substances 
should be considered as signaling agents involved in CDC clustering as well. We have 
investigated the role of acetaldehyde and ethanol |47j . and have shown that we can reset 
phase and induce autonomous oscillation with pulses of acetaldehyde. Further, we showed 
that both acetaldehyde and ethanol injections act differentially in different phases of the 
CDC and YAO and both act as agents of delay, supporting both our model and the role of 
these substances in processes. 

Finally, we note that our analysis is quite general and may find applications outside the 
realm of the yeast. We note that other micro-organisms can use signaling and coordination 
of various processes [H HH [I3l ETJ |35] . It is well recognized that feedback can produce 
synchrony in physical and biological systems and this has been extensively studied [H [T^ 
[Wl 1341 I45j (and many others). Clustering is a far less understood phenomenon. As early 
as 1977 a model of the CDC was introduced with periodic blocking at division [38], and 
clustering behaviour was also observed there. We note that clustering has been observed in 
Kuramoto type phase oscillators [43] and in certain all-to-all coupled networks of integrate 
and fire oscillators [58] and recently a rigorous proof of clustering and the stability of 
such solutions was obtained [36]. In that work the mechanism that leads to clustering has 
similarities with the idealized models in 
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A Appendix - Relation to PDE models 

Consider the hyperbolic conservation law: 

^ + ^{b{x,[u])u) = (A.l) 

for X on the unit circle, where u{x, t) > represents the density of cells of size x and [u] 
indicates that b depends functionally on the solution u i.e. b{x, [u]) depends on the values 
of u at all points x & at time t. This is a conservation law since f^i u{s, t) dx is a 
constant in time. This is the natural PDE version of our model of RS feedback (compare 
\12 \ \2Q \ [211 148j) and we provide it here for comparison with the existing literature. 
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If b has the form of RS feedback, then given any constant c > 0, there is a solution u{x) 
such that u{x) = c for x £ R'^ and u{x) for x G R is determined exphcitly by the feedback 
b{x, [c]), i.e., 

c, for X ^ R, 

for X G R. 



u{x) 



(A.2) 



Kb{x,[c]) 

This solution is analogous to the uniform solution in Proposition 2.3. 



Note that the corresponding version of (1.3) is given by: 

b{x, [u]) -- 

where 



1, if X ^ i? 
! + /(/), if xGi?, 



(A.3) 



I = u dx, 
Js 

i.e. the fraction of cells in S. For this feedback term, the uniform solution is piecewise 
constant with value c/(l + f{c\S\)) in R. 



We are not aware of any treatment of PDEs such as ( A.l ) with b of the form b{x, [u]). When 
there is no CDC feedback and diffusion is added to the PDE model (in various ways), then 
the uniform solution is asymptotically stable \12\ I20j. 



B Appendix - Analysis for k = 2 

Recall the notation a = /(^). 

We report the details of computations in the case where (1 + a)s < 1 — r. The other 
case can be treated similarly and we only give below the resulting expression of F. When 
(1 + a)s < 1 — r, there are 4 situations depending on the location of xi: 

• xi < r — s. In this case, xq leaves S before xi enters R. The point xi is not submitted 
to any feedback and hence xi{t) = xi + t for all t which implies F{xi) = ti = 1 — xi. 
(The occurrence of this case is independent ofr + (l + a)s<l.) 

• r — s < xi < r. Here xi is influenced, but not during the entire responsive region 
since xq gets out of S before xi reaches 1. More precisely, we hav^ 

{xi+t if < t < r - xi 

r + {1 + a){t — r + xi) if r — xi < t < s 
r + {1 + a){s — r + xi) + t — s if s<t 

It follows that -F(xi) = 1 — (1 + a)xi + a(r — s). 

• r<xi<l — (1 + a)s. Then 



xiit) 



xi + {l + a)t if 0<t<s 
xi + (1 + a)s + t — s if s<t 



*the occurrence of the late phase where xi, although being in R, moves with velocity 1, occurs due to 
the condition r + {1 + a)s < 1. 
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from which we obtain F(xi) = 1 — xi — as. 

• 1 — (1 + a)s < xi < 1. In this case, xi starts sufficiently close to 1 to have velocity 
1 + a when reaching the boundary. We have F{xi) = ^^(1 — xi). 

Recapitulating, in the case where r -\- (1 + a)s < 1, we obtain the following expression of a 
continuous decreasing map whose plot is given in Figure [s] (a) 



F{xi) 



1 — xi if < xi < r — s 

1 — (1 + a)xi + a{r — s) if r — s < xi < r 

1 — xi — as if r < xi < 1 — (1 + a)s 

^(1-xi) if l-(l + a)s<xi<l. 

In the case where r + (1 + a)s > 1, the map F is as follows (see Figure [s] (b)). 

1 — xi if < xi < r — s 



(B.l) 



F{x,] 



1 - (1 + q)xi + Q(r - s) if r - s < xi < + -^r - s 



2^1 + ^(1-0 if ^ + ^r-s<xi<r 
^(1 - xi) if r < xi < 1. 



(B.2) 



As argued for arbitrary k, the map F has a fc-periodic orbit which, for = 2, is composed 
of the boundaries and 1. By the Intermediate Value Theorem, it must also have a fixed 
point on the diagonal. 

The graph of F coincides with the anti-diagonal (1 — xi) for xi < r — s. For a = /(j) > 
it is strictly lower than this line if r — s < xi < 1; but if a < 0, then it is strictly greater 
than 1 — xi for r — s < xi < 1. The dynamics can be characterized completely for arbitrary 
parameters when k = 2. The following conclusions hold: 

• If r — s > ^, each point of the interval [1 — r + s,r — s] is part of a 2-periodic orbit 
X I— >• 1 — X I—)- X, except for the point x = ^ which is fixed. The return map F"^ thus 
has an interval of neutrally stable fixed points centered around xi = ^. 

• If Q > 0, if xi is above (resp. below) this neutral interval, we have F^(xi) > xi 
(resp. F^(xi) < xi) and so any initial point converges to 1 (resp. 0). 

• If a < 0, if xi is above (resp. below) this neutral interval, we have F^{xi) < xi 
(resp. F^(xi) > xi) and so any initial point converges to the boundary of the 
neutral interval. 

• If r — s = ^, there is a unique fixed point. It is stable for negative a and unstable for 
positive a. 

• If r — s < I there are three possibilities depending on where the diagonal line x = y 
intersects the graph of F (see Figure [s]). 

• If the diagonal intersects the second segment then there is a unique fixed point 
which is stable if a < and unstable if a is positive. 

• If the diagonal intersects the third segment of F, then there is again an interval 
of neutral period 2 points. The edge of the interval is stable for negative a and 
unstable for positive a. 
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Figure 8: Plots of the mapping F for k = 2 and a = /(^) = |. (a) r + |s < 1. (b) 
r+|s>l. Ifr — s>l/2 (which implies M > 2) then there is a neutral fixed point at 
f/2 that represents an isolated 2 cluster cyclic solution. If the second segment intersects 
the diagonal, then there is an isolated 2 cluster cyclic solutions that is stable for a < and 
unstable for a > 0. There also may exist neutral fixed points for r — s < 1/2 (M = 1) and 
certain conditions on the parameters where the third piecewise segment of F intersects the 
diagonal line. These fixed points represent 2 cluster cyclic solutions that are not isolated, 
but yet are neutral. 



• If the diagonal hits the boundary between segments 2 and 3 of F then there is 
a unique fixed point which stable for negative a and unstable for positive a. 

All of these possibilities are summarized in just four distinct types of behaviour in section 
6.2 and Figure [6j 

The condition r — s > 1/2 corresponds to M > 2. In such a case the cyclic 2 cluster solution 
consists of isolated clusters and it is contained in an interval of neutral period two points. 
This interval is an attractor for negative feedback and a repeller for positive feedback. 

Note that the condition r — s < 1/2 corresponds to \R\ + l^l > 1/2 which implies that 
M = [{\R\ + ISl)"^ \ = 1. Thus when the 2 clusters cannot be isolated, there is a cyclic 
2 cluster solution which is a fixed point of F. This solution may be unique and stable 
(negative feedback), unique and unstable (positive feedback), or neutral, depending on the 
parameters. 

One can easily calculate that the diagonal line cannot intersect the fourth segment of the 



graph of F in either case (B.l) or (B.2). 

In (B.l) it is seen that the third segment can be intersected by the diagonal by making r 
sufficiently small. This corresponds to the xi coordinate of the fixed point being greater 
than r so that the second cluster begins in the responsive region R. From the conditions. 
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it is still in R when xq leaves S. There is an interval of fixed points for F'^ even though the 
clusters are not isolated. 



The diagonal also can intersect the third segment for (B.2) if xq is in 5 when xi enters R, 
and when xi reaches 1, xq is still in S. Thus we have another case of interacting clusters 
that still leads to a neutral fixed point. 



C Appendix - Proof of Proposition 7.1 

Proof. The evolution of the system may be described qualitatively in terms of 'milestones,' 
e.g. Xq reaching s, or x^-i reaching 1. We first consider a system which evolves via the 
following sequence of events (which we call Case 1): 

Xk-i ^ r, xo i-> s, Xk-i ^ 1 (C.l) 

at which point after relabeling Xi := we have the initial condition again. By calculating 
the time taken in each step and finding the final value of xq, one can use the relation 
Xk-i = {k — l)d to solve analytically for d, getting: 

l + /3(r-.) 

Note that the Case 1 sequence will occur provided that ,s < xi = d and x^-i = {k — l)d < r. 



Using (C.2) in these relations gives (7.1). 



Two other sequences of events are possible and occur when we drop each of the constraints 



in (7.1). Case 2, when we allow r < (k — l)d, is characterized by: 

Xo ^ S, Xk-l H> 1, 
and Case 3, when we allow s > d, by: 

xi I—)- s, Xfc_i I— )• r, I—)- 1. 

Following the same procedure as in Case 1, we obtain for Case 2: 

l-sp 

and for Case 3: 

l + r/3 

d 



k{l + (3) 



A simple calculation shows that for a given M, cases 1, 2 and 3 exhaust the parameter set in 
(r, s) . We observe that there is never more than one cluster in the signaling region when the 
response region is nonempty, and thus the dynamics of the system are determined entirely 
by/5 = /(|)- A graphical representation of the regions of parameter space corresponding 
to these three cases can be seen in figure ([7|). 

The map F is affine in a neighborhood of the fixed point, i.e. x i— )• Ax + b where x = 
(xo, ...,Xfc_i)"^ and A is a matrix. We next analyze A in the three cases. 
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Case 1: 



A 





1 
1 



-(1 + /?) 



••• + 

••• 1 -{1 + 13) 

Thus we may determine the stabihty of the fixed point by studying the eigenvalues A^ which 
has characteristic equation 



A 



fc-i 



(1 + /3)(A'=-2 + A'^-^ + --- + A + 1 



fc-3 



0. 



(C.3) 



Notice that A = 1 can easily be ruled out as a root. For A 7^ 1 we can rewrite (C.3) as 

n-l 



' i=0 



A" A" - 1 

+ - = 0. 



1 + 



A-1 



After simplification, we see that A 7^ 1 is a solution of (C.3) if and only if 

■A + /3^ 



1 + 



A" = 1. 



(C.4) 



Now, suppose ^ > 0. If |A| < 1 (and A / 1), then |A"| < 1 and |A + /3| < |1 + Thus 

A + /3 



1 + 13 



|A"| < 1, 



i.e. A cannot satisfy (C.4). Thus for positive feedback all of the eigenvalues lie outside the 
unit disk and so the map is unstable at the fixed point. Further, it is not only unstable, 
but is unstable w.r.t. all possible perturbation directions. 

For the case /3 < 0, suppose | A| > 1. Write A + /3 as A — (— Then by the reverse triangle 
inequality, |A + /3| = |A-(-/3)| > ||A| - | - /3|| = ||A| = |A| + /3 > 1 + /3. Thus we have 

A + /3 



1 + 13 



|A"| > 1, 



and ( |(I4| ) is not satisfied. Also, if |A| = 1 but A / 1, then |A + /3| < 1 + 13. Therefore 
if /3 < then all the eigenvalues of A lie on the interior of the unit disc, and the map is 
stable. □ 

The stability results for case 1 are illustrated in Figure|9]for k = 2,4, ... ,12 and /? = /(^) G 
(-.5, .5). 

Case 2 and Case 3: 

In cases 2 and 3, the linear part of the map at the fixed point is represented by the matrix: 



A 





1 
1 












1 



29 



Spectral radius / minimum eigenvalues 



1.10 

1.05 

1.00 

0.95 

0.90 
1.10 

1.05 

1.00 



0.95 - 



0.90 





/ 


,/ 




k = 2 




0.0 0.2 0.4 










k = 8 





1.10 

1.05 

1.00 

0.95 

0.90 
1.10 

1.05 



1.00 



0.0 0.2 0.4 



- 0.95 



0.90 











k = 4 




0.0 0.2 0.4 










k = 10 





1.10 

1.05 

1.00 

0.95 

0.90 
1.10 

1.05 



1.00 



0.0 0.2 0.4 



- 0.95 



0.90 















k = 6 




0.0 0.2 0.4 










k = 12 





0.0 0.2 0.4 



Figure 9: Case 1. The spectral radius / smallest eigenvalue modulus for k = 2,4:, ... ,12. 
The X-axis is the feedback parameter /3. For negative /5, the y-axis is the spectral radius of 
A. For positive (3, the smallest eigenvalue (w.r.t. modulus) is plotted. The plots show that 
the k = M + 1 solution is stable for negative feedback and unstable (in all directions) for 
positive feedback. Notice that as k grows, the stability/instability becomes weaker. 

This matrix has characteristic equation A'^^^ + A'^'^^ + ... + A + 1 = whose roots all have 
absolute value 1. Thus the map is neutrally stable in both cases 2 and 3. □ 



30 



